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Abstract 

This article discusses a potential impact of turbulent velocity fluctuations of the air on a drizzle formation in Cumulus 
clouds. Two different representations of turbulent velocity fluctuations for a microphysics formulated in a Lagrangian 
framework are discussed - random walk model and the interpolation, and its effect on microphysical properties of the 
cloud investigated. 

Turbulent velocity fluctuations significantly enhances velocity differences between colliding droplets, especially 
those having small sizes. As a result drizzle forms faster in simulations including a representation of turbulence. 
Both representations of turbulent velocity fluctuations, random walk and interpolation, have similar effect on droplet 
spectrum evolution, but interpolation of the velocity does account for a possible anisotropy in the air velocity. 

All discussed simulations show relatively large standard deviation of the cloud droplet distribution from the 

onset of cloud formation is observed. Because coalesence processes aerosol inside cloud droplets, detail information 
about aerosol is available. Results from numerical simulations show that changes in aerosol spectrum due to aerosol 
processing during droplet coalescence are relatively small during ~2Q min. of the cloud evolution simulated with 
numerical model. 

Drizzle forms initially near the cloud edge, either near the cloud top, where the mass of water is the largest, or near 
the entrainment eddies. 
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1. Introduction 


One of the unresolved problems in cloud physics 
is drizzle formation in a warm (ice free) clouds. 
Ear ly mode ll ing st u dies using Lagrang i an par cel mod¬ 


els Warner ( 1969h . Bartlett and Jonas (1972) demon¬ 


strated that condensational growth leads to a very nar¬ 
row droplet spectrum, contrary to the observations, where 
the d roplet spectrum inside the cloud was relatively broad 
e.g. IWarnei ( 1969 ). A narrow droplet spectrum makes 


'Now at tlie: Met Office, FitzRoy Road, Exeter, EXl 3PB, United 
Kingdom, miroslaw.andrejczuk@metoffice.gov.uk 


collision between droplets inefficient and as a result it 
takes a long time to switch from a condensational to 
a coagulational droplet growths. Observations show, 
that develop ment of pre c ipitati on in the clouds may be 
rapid e.g. Poke et all (2007) which couldn’t be ex¬ 
plained by the modelling results. Several mechanisms 
to explain this discrepancy between numerical model 
and observations were proposed and were reviewed by 
Beard and Ochsl l 1993h . 


Undoubtedly representation of the cloud formation pro¬ 
cess by a parcel model is an approximation. Nevertheless 
attempts were made over the time to explain the width of 
the cloud droplet spectrum in clouds using this approach. 
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Yum and HudsonI (120051) showed that a much broader 
droplet spectrum can be obtained when instead of one ver¬ 
tical velocity, a velocit y distribution is used. Another ap¬ 
pro ach was utilized by IBewlev and Lasher-Trapnl (120051) 
and Lasher-Trapp et al. ( 20051) . who used a velocity field 
from a bulk model to derive trajectories for the par¬ 
cel models. For both these approaches droplet spectra 
broader than predicted by a single parcel were reported, 
however, in these approaches the microphysics is sepa¬ 
rated from dynamics and thermodynamics of the Eulerian 
model. 

The use of th e bin models iGrabowski 

(Il99lh 


Kogan 


Ackerman et al 


(fi^. 


Feingold et al 


Khain et al 


LeroY.e^^ (120091). lOvchinnikov and Easteri (120101) . 


Blvth et al.l (120131) . IWvszogrodzki et al.l (120131) . where 


1989) 


1994), 


(2004), 


droplet spectrum is represented as a continuous function, 
made the problem of the transition from condensational 
to coagulational disappear. Even very high resolution 
(in bin sizes) parcel models with bin microphy s ics ar e 
capable of producing pr ecipitation ICooper et al.l (Il997h . 
Pinskv and Khain ( 20021) independent of whether turbu¬ 
lent enhancement is or is not taken into account. It’s not 
clear, however, whether the rain forms in these models 
due to physical processes or numerical errors associated 
with the numerical solution of the condensational growth 
and collision or the physics itself, since detailed com¬ 
parisons between Eulerian (bin) and Lagrangian models 
were never to our knowledge reported. 

In the recent years a new approach to microphysics 
formulated in Lagra ngian framework was p r opose d 


for both warm ra in lAndrejczukeLal 


Shima et'd] ( 200 91). Riechelmann et al. _ _ ___ 

clouds Solch and Karcher ( 201ol) . In this approach dy¬ 



grangian parcels overcomes many problems present in 
Eulerian models. The Lagrangian representation of mi¬ 
crophysics is diffusion free; each parcel can be treated 
individually, which makes representation of the sub-grid 
scale variability easier; the edge of the cloud is re¬ 
solved without the need of the us e of a special technique s 
(for instance VO F discussed in iMargolin et alJ (Il995h . 


Kao et alJ (119991) ). At the same time representation of the 


field with a limited number of parcels may lead to random 
fluctuations in derived fields (for in stance concentratio n 
of droplets and aeros ol cloud water) Kenzelbach ( 1990l) . 
Salamon et al.l ( 20061) : which in an extreme cases may 


lead to parcel free computational grids. 

This article focuses on the representation of the turbu¬ 
lence in a numerical model with a Lagrangian represen¬ 
tation of microphysics and its effect on drizzle formation. 
The effect of air turbulence on drizzle formation is a com¬ 
plex and not fully understood problem in cloud physics. 
Turbulence can affect directly relative velocities of the 
colliding droplets, leading to a coalescence of droplets 
which would not collide in a laminar flow. It can also indi¬ 
rectly influence drizzle formation, by modification of the 
environment in which droplets grow/evaporate. A broad 
review of the effect of tu rbulence on cloud s and its impor ¬ 


tance were discuss e d by: Srivasta^ ^^98^ .Tonasid 9961) 
Bodenschatz et alJ (l201()l). Devenish et al. ( 2012 ) 


and 


Grabowski and Wangl (12013 ). 

The effect of air turbulent velocity fluctuations 
droplet motion in bin models may be rep¬ 
an enhancement of collision efficien- 
a high resolution turbulence sim- 


on 

resented as 
cies derived from 
ulations 


Franklin et al. (2005), Pinskv et^ ( 20061) . 


Wang et al. j 20061). This m e thod has been used recentl y 


namics and thermodynamics are represented in the tradi¬ 
tional Eulerian framework, whilst microphysics is repre¬ 
sented in a Lagrangian framework with two way interac¬ 
tions between Eulerian and Lagrangian parts. Lagrangian 
microphysics tracks La grangian parcels (s ometimes re¬ 
ferred as super-droplets Shima et al. ( 2009h ). each repre¬ 
senting a number of real aerosol, having the same chem¬ 
ical and physical properties. Depending on the condi¬ 
tions determined by the Eulerian model, water can con¬ 
dense on or evaporate from the aerosol surface. Result¬ 
ing forces, together with a drag force are return to the 
Eulerian model. Transport of physical properties by La- 


bv lPinskv and KhainI (l2002h . IWvszogrodzki et al.l (l2013l) . 
Also, a recently reported model with Lagrangian micro- 
phvsics lRiechelmann et al.l (l2Q12l) uses this approach. Be¬ 
cause in the Lagrangian microphysics the parcel velocity 
is a predicted variable, it can be used with the parameteri¬ 
zation of the sub-grid scale velocity fluctuations to deter¬ 
mine the relative velocity of the colliding parcels. 

This article reports an application of a two possible rep¬ 
resentations (parameterizations) of turbulence for a model 
with microphysics formulated in a Lagrangian frame¬ 
work, one - random walk model, and the other interpo¬ 
lation of the velocity to parcel location and investigates 
its effect on a droplet spectrum. The results from these 
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models are compared with the results from the model 
where parcel terminal velocity is used when calculating 
probability of droplet collisions. When for a parcel pre¬ 
dicted by the model velocity is used to calculate proba¬ 
bility of droplet collisions it is assumed that droplet tur¬ 
bulent transport and interaction occur in the same scales, 
because the same velocity is used for both. Assump¬ 
tion that colliding droplets have terminal velocity corre¬ 
sponds to the case, where turbulent transport and interac¬ 
tion (collision-coalescence) between droplets happen in a 
different scales. This gives a lower and an upper limit of 
the effect of the turbulence on droplet collisions, because 
in former case it is consistent with the turbulent transport 
model and in latter is neglected. 

The next section describes numerical model and repre¬ 
sentation of turbulence. Results are discussed in section 
3, and conclusion are in the last section, 4. 


2. Numerical Model 


2.1. Eulerian Model 


In this article the 
as a driving model 
physics. Eulag 

solver ([Smo larkiewicz and Margolin 


EULAG model is used 
for Lagrangiam micro- 
an Eulerian/Sem i-Lag rangia n 


Prusaetal 


(!2008h ). with the Eulerian 
the simulations reported in this article. 


(Il997h . 
version used 
Model 


-I- Vpviff = + pF^ + pFp (1) 


in tne simulations reportea in 
equations in an anelastic approximation can be written 
as: 

dpiff 
dt 

where p is density, diffusion coefficient, ij/ - any de¬ 
pendent variable (u, w, 0, qv) with its associated Eule¬ 
rian forcing - F^, Fp describes forces from a Lagrangian 
model and are dehne d below. Km is de rived from a prog¬ 


nostic TKE equati on iDeardorfli (Il973h. see details about 
implementation in Margolin et al.l 'i 1999 ). 


2.2. Fagrangian Microphysics 

Lagrangian microphysics tracks Lagrangian parcels, 
each representing a group of real aerosol particles. Each 
parcel can be characterized by the same velocity and size 
(both dry aerosol radius and radius of the droplet if water 
condenses on the parcel), and occupy approximately the 
same place in a physical space. 


Lagrangian model equations describing evolution in 
time of the droplet radius (r,), velocity (v,) and location 
(xi) for a parcel p can be written as; 


dr, 


du 


-p- = + Up.) + 6i2g 

dt Ti ‘ 

dxp 


^ = -(5 

dt rp 


1 


■ S eq) 


( 2 ) 


(3) 


(4) 


where rp -droplet radius, S - supersaturation at parcel lo¬ 
cation, S eq- equilibrium supersaturation for given aerosol 
size and temperature (see lAndrejczuk et al.l(l2008h for de¬ 
tail description of other symbols). The u‘^ - deterministic 
air velocity and u - turbulent component determined from 
the sub-grid scale mo del; g is gravi t y. Equ ationlUis solved 


using VODE solver (IBrown et al.l 119891)), equation [3] i 


solved using the backward Euler method, and equation |4] 
is solved using the forward Euler method. 

One of the established ways to represent turbulent com¬ 
ponent of the flow velocity is by a random process (ran¬ 
dom walk) having normal distribution and standard devi¬ 
ation yjlKi!At, which is derived in the Appendix. Note 
that to determine Ki different mixing length than to deter¬ 
mine K„, in eq. [T]was used. Although the mixing length 
of the order of model grid size is typically used in sub¬ 
grid scale model, in the Lagrangian model because the 
location of each parcel within the grid is known, the mix¬ 
ing length based on TKE and dt can be derived. This, 
simplihed, description of the turbulence treats turbulence 
as a random walk process and represents the effect of air 
turbulent velocity on Lagrangian parcels velocity through 
the variability in the Eulerian model velocity held within 
a computational grid. A broad review of theory and ap- 
plication of Lagrang i an tra nsport models can be found in 
Wilson and Sawford ( 19961) . No attempt has been made to 
represent eddies present in the turbulent how, but not re¬ 
solved by the model; it is assumed that at each time-step 
turbulent Huctuations of the air velocity are independent. 
It is achieved by generating a new random number for 
each Lagrangian parcel on each time-step. This me thod 
follows an approach proposed by De Almeida ( 1976h . 

The other method of representing variability of the ve¬ 
locity within the computational grid is to use interpola- 
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tion; and instead of using the same value of the velocity 
for all parcels within a given computational grid the ve¬ 
locity can be interpolated from neighboring grids to each 
parcel location. 

These two approaches, later referenced as a turbulence, 
provide a way to describe the variability of the air turbu¬ 
lent veloc ity inside comput ational grids in the spirit dis¬ 
cussed by Srivastava ( 19891) for a microscopic supersatu¬ 
ration. 

Details of representation of the coalescence process 
in Lagrangian m i croph ysics is described in detail i n 


Andreiczuk et al. (l2O10l) and Andreiczuk et al. ( 2012 ). 


Collisions of all Lagrangian parcels within the collision 
grid (this does not haye to be the same as Eulerian grid, 
in the simulations reported in this article collision grid 
is a quarter of a computational grid), with sizes larger 
than 3 ti m and haying water on the surface are consid¬ 
ered, and iLond (119741) analytical expression for grayita- 
tional collision efficiency is used. Each collision eyent, 
based on the size of colliding droplets and aerosol size in¬ 
side the droplets is assigned to one of the pre-defined mi¬ 
crophysical grids spanning both aerosol sizes and droplet 
sizes. Eor each microphysical grid the mass of the aerosol, 
the mass of the water and the number of new droplets 
is calculated. Based on this information new parcels are 
created for each microphysical grid for which the num¬ 
ber of physical droplets is larger than specified thresh¬ 
old leyel - 156.25 in simulations discussed in this arti¬ 
cle (which corresponds to resolying 1 droplet/m^). Newly 
created parcels for collision grids at the edge of the cloud 
are placed next to randomly chosen parcels within this 
grid. Eor grids within the cloud they are randomly placed 
within the Eulerian computational grid. Additionally the 
algorithm assures also that existing parcels represent a 
larger than threshold leyel number of physical particles. 
Should the collision lead to a smaller number, the prob¬ 
ability of collision, with all other parcels for this partic¬ 
ular one is reduced. This representation of coalescence 
process allows to process not only droplet sizes, but also 
aerosol sizes during droplet coalescence. It is assumed 
that newly created Lagrangian parcels moye with termi¬ 
nal yelocity. 


2.3. Models coupling 

Drag force and tendencies for a temperature and water 
yapour mixing ratio are calculated in Lagrangian model 


and used in an Eulerian part; 


P e g>'>du.j.k) ^ 

(5) 

P e gndpj.k) 

(6) 

Oe L 

Fpe = Fpq„ 

^ e ^ p 

(7) 

where indexes i, j, k index numerical model computa¬ 
tional grids, Mj is the number of real particles parcel i 
represents. A Tg and 6^ are temperature and potential tem¬ 
perature profiles. Treatment of these forces is similar to 
the treatment of the sub-grid scale tendencies as discussed 

by Margolin et al. ( 19991). A finite-difference ai 

^proxima- 

tion to eq. fTIcan be written as (Margolin et al. (1 

19990: 

iAr‘ ^ ADV+0.5AtFf\ 

(8) 


where ij/ - t)/" + 0.5AtF” + AtD^ + AtFp^, are 
forces associated with pressure gradient, a bs orbers , 
and buoyancy ( Smolarkiewicz and Margolinl ( 1997 ). 


Margolin et al. ^1999l) l. 


yection, 

DATA 


which 

scheme 


ADV 

calculated 


denotes 

using 


ad- 

MP- 


_ _ _ _ (ISmolarkiewicz 

Smolarkiewicz and Grabowskil ( 199(ilK 


(119831) . 


2.4. Model Setup 


In this article a 2D idealized setup is used to inyes- 
tigate the effect of air turbulent yelocity on drizzle for¬ 
mation. Initially, an atmosphere with a IjOdOldz - 1.3 
nT^ was specified. Below 2km a relatiye humidity was 
defined as 85%, dropping to 75% aboye this leyel. The 
model domain coyers 3.2 km x 5.0 km, resolyed with 
25m resolution in each direction. Three-modal, log¬ 
normal aerosol distribution, corresponding to continental 
air, h aye been specifi ed for the whole domain using yalues 
from Whitby ( 1978 ): r„ =0.008 pm, cr„=L6, A„=1000, 
ra =0.034 pm, cra=2.1, Na- 800, rc =0.46 pm, crc-2.2, 
Nc^O.72. 

An initial aerosol distribution in each computational 
grid was represented by 100 parcels. The aerosol spec¬ 
trum used in the simulations is shown in figure [T] Aside 
from a spectrum ayeraged oyer the whole domain, the 
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standard deviation of the values for each bin is also plot¬ 
ted. Variability in each bin is a result of a random 
sampling of the initial distribution with finite number of 
parcels. Although on average the spectrum agrees with an 
analytical distribution, sampling with a limited number of 
parcels leads to a variability of the distribution for a dif¬ 
ferent grids and as a result it affects droplet number and 
size in each model grid during the condensational growth. 

The model was forced by the surface temperature 
source prescribed as: 


f,(x, z) = Athcxpi- 


(x - xq)^ + {z- Zq)^ 
800 


), (9) 


with A,/,=0.15 K/s and (xq, zq) set to (0,400). Simula¬ 
tions have been run for 1080s, with a time-step 0.25s. Mi¬ 
crophysical grid in aerosol space (in pm) was specified 
as: Kaik) - iQ-^- 2 +k*o. 2 ^ radius space 28 

bins were used, and radius for a bin k was specified as: 
r(k) - p^l^rik - 1), with r(l)=ljum and p-2. 

The following model setups are discussed in details in 
this article: 


REF - reference simulation. No representation of tur¬ 
bulence in eq. [3(v'=0). Deterministic flow velocity 
u'^ at parcel location is determined by interpolating 
velocity from 4 (8 in 3D) nearest Eulerian computa¬ 
tional grids to a parcel location. The vertical compo¬ 
nent of the parcel velocity is used to evaluate colli¬ 
sion kernel. 


BIN - As far as parcel movement is concerned it is as¬ 
sumed that parcel velocity is equal to the air velocity 
interpolated to a parcel location, that is eq. |3]is not 
solved. The terminal velocity of the parcel is added 
to the vertical velocity. The parcel’s terminal veloc¬ 
ity ha^_been_toenn^d from the expression given 
by Simmel et al.l ( 2002 ). To simplify calculations in¬ 
stead of the terminal velocity for a parcel size, ter¬ 
minal velocity for the center of the collision bin to 
which the parcel is assigned is used. When the col¬ 
lision kernel between parcels is calculated it is as¬ 
sumed that the parcel velocity is equal only to a ter¬ 
minal velocity. This setup mimics bin model as far 
as collision-coalescence process is concerned and as¬ 
sumes that turbulent transport and collision happen 
in different scales (are independent). 


TURB - for a given parcel deterministic flow veloc¬ 
ity has 2 components deterministic: u‘^ - u(x) - 
^■^(pKi), where u(x) is Eulerian velocity for the 
computational grid to which parcel belongs, K/ - dif¬ 
fusion coefficient for a Lagrangian model, and p - 

gas density; and turbulent component v'= with 
the derivation of of this model presented in the Ap¬ 
pendix. 

REF2 - Similar to REF, but in this case in additional to 
the difference in vertical velocity also difference in 
horizontal velocity is taken into account when eval¬ 
uating the collision kernel. 

3. Results 


3.1. Evolution in space 

Eigure EJi-E};, shows snapshots of a cloud water 
mixing ratio {pc) for the model solution for times 8 
min. (fig. |^), 15 min (fig. |2]r) and 17 min (fig. 
|2t) for the REF run. Cloud evolution exhibits fea¬ 
tures typically found in a 2D Cumulus developing in 
a stable stratified atmosp h ere and reported alre a dy by : 


)wskil ( 


Klaassen and Clar kI (119851). iGrabowski and ClarkI (119911) . 


Brenguier and Grabowskil(ll993h . After reaching conden 


sation level water vapor condenses and forms the cloud; 
air continues moving upward due to the temperature ex¬ 
cess within the thermal. During the ascent in a sta¬ 
bly stratified air cloud mixes with the environmental air, 
forming entrainment eddies. The cloud water mixing ra¬ 
tio reaches values as large as 12 g/m^ during the simula¬ 
tion, and at the end of the simulation the maximum ver¬ 
tical velocity is around 18 m/s. The cloud water mixing 
ratio shows relatively large variability in space, being a 
result of the fluctuations of the parcel number within the 
computational grid and also because of the variability in 
the randomly generated aerosol spectrum for each model 
grid. 

Because in the numerical model the coalescence pro¬ 
cess is present, coalescence of the droplets in time leads 
to a drizzle formation. Eigure |2}i and|2^ show evolution 
in time of the q^ - the rain water mixing ratio for a REE 
case. A 50 pm droplet radius is defined to be a border 
between a cloud and a rain droplet sizes. Large droplets 
reside near the edge of the cloud, and at the end of the 
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simulation qr exceeds 3 g/kg. This behavior is typical for 
the three setups; REF, TURB and REF2. For a BIN case 
negligible forms. 

The edge of the cloud is a place where the difference in 
the vertical velocity within the collision grid is the largest 
(either because of the gradient of velocity between the 
interior and exterior of the cloud for REF or REF2, or 
because TKE, used to derive stochastic velocity pertur¬ 
bations is largest). Additionally at the edge of the cloud 
the droplet spectrum can be broader than in the core of the 
cloud, which also enhances coalescence between droplets. 
The width of the droplet spectrum increases near the edge 
of the cloud because: a) the gradient of the supersatura¬ 
tion is largest there and as a result the interpolation of the 
thermodynamical parameters to a parcel location within 
the same computational grid different droplet populations 
may grow and evaporate at the same time depending on 
the distance from the edge of the cloud; b) vertical veloc¬ 
ity changes sign near the cloud edge (in figurel^the solid 
line shows a contour of the value 0), as a result droplets 
evaporate in a down-drafts near the cloud edge; and c) 
because largest droplets, formed in the center of the up¬ 
draft, after reaching cloud top start moving along cloud 
edge. 

3.2. Evolution in time 

Evolution of the TWP (Total Water Path) and RWP 
(Rain Water Path) in time for cases discussed are shown 
in figure [3] During the early stage of the cloud develop¬ 
ment TWP is very similar for all cases. In time, how¬ 
ever, differences appear as a result of the interaction of 
the cloud and flow dynamics. The largest amount of wa¬ 
ter is for REF2 and REF and smallest for TURB. The 
decrease in a TWP for a TURB case can be associated 
with a random velocity of the air within the grid. As a 
result evaporation at the edge of the cloud may be larger 
for this case, because droplet trajectories can be different 
than for other three cases when velocity is interpolated to 
a parcel location. Much larger differences are observed 
for the RWP. The largest amount of water is for the REF2 
case, the simulation where in coalescence calculation full 
velocity is taken into account. Compared to REF simula¬ 
tion, REF2 has around 15% more RWP. Eor the BIN sim¬ 
ulation negligible amount of RWP has formed during the 
simulation time. But undoubtedly the coalescence process 
is active for this case also and droplets as large as 40 


have been formed in this simulation. These droplets are 
much smaller than for the three other simulations, where 
sizes up to 370 are present. 

3.3. Vertical profiles 

Vertical profiles of the mean radius (r), standard devia¬ 
tion of cloud droplet distribution (cr), cloud water mix¬ 
ing ratio (qc) and number of cloud droplets are shown 
in figures |4] for a REF case (with very similar statistics 
observed for a REF2 and TURB simulations) and for a 
BIN case in figure |5]for a 3 times: 6, 8, and 15 minutes 
from the beginning of the simulation. These 3 times show 
cloud characteristics at the very early stage of the cloud 
development (the time when qc becomes larger than 0), 
through the initial formation of the drizzle (space distri¬ 
bution after 8 min. is shown in figures |2]r and|2}l), to the 
moment when significant drizzle develops (plots |2j; and 
|2^). These diagnostics were calculated for each model 
level by taking into account only model grids, where qc 
was larger than 10“^ g/kg. Eigures |4] and |5] show a very 
similar development of the cloud for both cases. Initially 
a small amount of water condenses, but at the same time 
the mean radius is around 5 pm and standard deviation 
around 1 pm. Around 1/3 of the total aerosol concentra¬ 
tion activates at this stage. In time, the cloud thickens and 
after 8 minutes the cloud base moves upward. The mean 
droplet radius increases with height and reaches 12-13 pm 
near the cloud top after 15 minutes At that time the num¬ 
ber of cloud droplets decreases with height for a REF case 
from 500 cm“^ to 250 cm“^. Eor the BIN simulation, the 
cloud droplet concentration does not change significantly 
with height and oscillates around 500 cm“^. Another dif¬ 
ference between these 2 cases is in the standard deviation 
of the cloud droplet distribution after 15 min. above 2.5 
km, and for the BIN case standard deviation is around 1.5 
pm, except at the cloud top, where it reaches 4.2 pm. Eor 
REF case it is of the order of 3.5 pm, with the value near 
the cloud top of ~5 pm. At earlier times and below 2.5 
km standard deviation is very similar for both cases. 

3.4. Droplet and aerosol spectra 

Because for each parcel information about droplet size 
and aerosol size is available, the relation between aerosol 
sizes and droplet sizes can be derived. Eigure|6]shows this 
relation mapped on a Eulerian microphysical grid after 18 


6 


minutes. The relation between aerosol size and droplet 
size is complex and for a given aerosol size there is a 
broad range of droplet sizes formed on it. Because ini¬ 
tially aerosol sizes were limited to 1 fim, sizes larger than 
that were created in a coalescence process. Much larger 
sizes of both aerosol and droplet are present for a REF 
case, where the coalescence is more intensive. 


From observations or spectrum resolving (bin) models 
typically information about either aerosol distribution or 
droplet distribution is derived. This information for mi¬ 
crophysics in a Lagrangian framework is obtained by in¬ 
tegrating relation shown in figure |6] along one of the di¬ 
mensions. Averaged over the whole cloud the droplet 
spectrum is shown in hgure |2l After 9 minutes - figure 
|7^ the droplet spectrum is very similar for all cases. In 
time, similar to other cloud properties the differences be¬ 
tween BIN and REF, TURB and REF2 emerge (|7}r,|7]:, 
EH). Formation of the large droplets for the BIN simula¬ 
tion is much slower than for the remaining three cases. In 
an early stage of cloud development the largest droplets 
form for the TURB case. At the end of the simulation, the 
spectra fo^:REF, TURB wdREFl cases are very similar. 


Evolution of an aerosol spectra averaged over the whole 
cloud is presented in hgure[8] Because coalescence is pro¬ 
cessing not only droplet sizes, but also aerosol sizes, in 
time the aerosol spectrum also changes. There are indica¬ 
tions, especially for later times, that aerosol is processed 
during droplet coalescence, but the differences between 
BIN (where aerosol processing is negligible) and other 
three simulations are small. Fastest aerosol processing is 
for a TURB simulations, with a smaller, for a REF and 
REF2. Much more time is needed for the coalescence 
to process aerosol and form large/giant aerosol than the 
length of an idealized simulation discussed in this arti¬ 
cle. The re is, however, evid ence that this aerosol process¬ 
ing (eg. ITwohv et al.l (l2013h ') may be an important source 
of large aerosol, because within the ~2Q minutes max¬ 
imum aerosol size has doubled, reaching 2.5 /urn radius 
(see hg. |6H) compared to 1 fxm initially. The increase of 
the concentration of the large aerosol is at the expense of 
the aerosol having sizes between 0.01 and 0.2 jum. 


3.5. Velocity enhancement 

The collision kernel, describing the probability of the 
coalescence of two colliding droplets can be written as; 

Kin, rj, Vi, Vj) = Mn + rjf\vi - Vj\Eij, (10) 


where r,.,rj is the radius of a droplet //j with the corre¬ 
sponding velocity v,7vj, Eij is the gravitational collision 
efficiency. With the equation of motion solved for each 
parcel in the model, diagnostics of the relative velocity of 
colliding droplets can be determined. Equation[T0]can be 
rewritten in the following way: 


Kin, rj, Vi, vj) = n{ri + rjY 


\Wi - Wj 


■ W/IE; ;, 


( 11 ) 


with the w,/wy being the terminal velocity for a parcel // j. 
An enhancement of the gravitational collision efhciency 
due to the turbulent velocity fluctuations can be dehned 
as: E,, =< >b, where vjvj is the velocity of col¬ 

liding parcels, either vertical only or vertical and horizon¬ 
tal for a REF2, and <>t - averaging operator. Although 
it is possible to calculate this enhancement for each pair 
of droplets it is easier to use a bin structure the same as 
was used to map collisions between parcels, and in such 
a case w,/wj represent a terminal velocity for the center 
of the bin rather than for an individual parcel, and aver¬ 
aging is done for each bin. Eigure |9] shows a enhance¬ 
ment of the gravitational collision efficiency. Eour panels 
show Ej, for a different bin sizes: 2.6 fim (microphysical 
bin 4) - figure |9H, 9.1 ^m (microphysical bin 10) - hg¬ 
ure |9H, 28.8/im (microphysical bin 16) - hgure |9]:, 94.3 
jum (microphysical bin 22) - hgure |9}l. The largest en¬ 
hancement is for small bins with a value reaching 125 for 
a TURB simulation. Simulation REF2 has value 120 and 
REE almost 105. Expectedly, the BIN simulation has a 
constant value of Ev - 1 independent on size, because the 
vertical velocity in this case is a droplet terminal veloc¬ 
ity. With the increasing droplets sizes the departures of 
enhancement factor from 1 is smaller, and for bin 22 fig- 
urefTOH it approaches 2. Eor a given bin, the largest E^, 
is found for adjacent bins. Eor bins separated by a large 
distance E,, is much smaller. The distribution of the E,, 
is non-symmetric, with larger values of E,, found for the 
sizes smaller than bin under consideration. The results 
presented in hgure |9] show that even small differences in 
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air velocity can affect velocity statistics for droplets hav¬ 
ing similar sizes and as a result moving with a similar ve¬ 
locities. For droplets having large differences in sizes tur¬ 
bulent velocity would have to be significantly larger than 
terminal velocity of larger droplets to influence E^,, which 
is not a case in the simulations discussed in this article. 
The large values of the enhancement of the gravitational 
collision efficiency for small droplets are assiciated with 
the fact that those droplets adjust quickly to the environ¬ 
mental flow, and as a result the velocity of small droplets 
is very similar to the velocity of the air. Much smaller val¬ 
ues of the enhancement of the gravitational collision effi¬ 
ciency for large droplets are because these droplets need 
more time to adjust to the environmental flow and their 
terminal velocity is of the same order as air velocity fluc¬ 
tuations. 

Diagnostics show that the maximum turbulent velocity 
for a TURB case is around 3 m/s and has the same order 
as a terminal velocity of the 650 fim droplet —2.5 m/s. 
Note, however, that 3 m/s represents the largest value and 
standard deviation of the turbulent velocity for the whole 
cloud is ~0.1-0.25 m/s. An example of the velocity fluc¬ 
tuations for run REF and TURB are shown in figure [TO] 
There is a significant difference between the statistics for 
these two cases. For a TURB case distribution of u and 
w are almost identical. For the REF case anisotropy of 
the velocity fluctuations is observed, and tails for a w dis¬ 
tribution is much broader than for a u , with the statistics 
for case laying in-between M and vr for a REF. 

4. Conclusions 

In this article two possible representation of the air tur¬ 
bulent velocity in a numerical model with a Lagrangian 
representation of microphysics are discussed. Air turbu¬ 
lent velocity in the model is represented either as a ran¬ 
dom walk process, with the standard deviation of veloc¬ 
ity fluctuations derived from the diffusion coefficient pre¬ 
dicted by the numerical model (Eulerian part) or as an 
interpolation of an air velocity to a parcel location. The 
random walk model is derived for an anelastic approxima¬ 
tion and it’s shown that an additional, deterministic, term 
IjpVipKm) needs to be included in a random walk model 
for consistency with the Eulerian model. It is argued that 
the mixing length used in a Lagrangian model should be 
smaller than used to derived diffusion coefficient for an 


Eulerian model; and the mixing length scale based on 
TKE and model time-step for use in Lagrangian micro¬ 
physics is introduced. It is demonstrated that interpolation 
of the velocity to parcel location and use of this velocity 
in a collision kernel has a similar effect to representation 
of the turbulence as a random walk, and can be treated as 
an alternative to a random walk model. Additionally, for 
turbulence treated as a random walk, because all parcels 
within the particular computational grid have the same 
value of the deterministic velocity, fluctuations in the par¬ 
cel number in a computational grid may be larger than 
for the case when velocity is interpolated to a parcel loca¬ 
tion. Unlike the random walk model, interpolation takes 
into account possible anisotrop y in the flow veloc i ty, als o 
observed in a laborat ory studies Malinowski et al.l ( 2008 ). 


Korczvk et al.l (120121) but in a much smaller scale. 


If the sub-grid scale transport and droplet collisions 
occur on the same scale, air turbulent velocity can sig¬ 
nificantly enhance velocity differences between droplets, 
especially those having small sizes, when the turbulent 
velocity is much larger than terminal velocity of these 
droplets. In the cases discussed in this article enhance¬ 
ment as large as 120 has been observed for the case 
when only vertical velocity is taken into account when 
calculating relative velocity of colliding parcels. Allow¬ 
ing for the differences in a horizontal velocity in¬ 
creases by ~15 %. The turbulent enhancement of the 
velocity of the colliding droplets obtained in this article 
is much larger than obtai ned form a DNS simulatio ns 
and used in EES mode l bv lWvszogrodzki et alJ (120 13h or 
iBenmoshe et alJ ( 2012 ). This difference is associated with 
the fact that in Lagrangian microphysics the same veloc¬ 
ity is used for the transport and when evaluating collision 
kernel. The values obtained with air turbulent velocity 
fluctuations {REF, TURB) provide an upper limit of the 
impact of the air velocity fluctuations on gravitational col¬ 
lision efficiency, with the lower limit given by values from 
the BIN simulation. An additional parameterization is 
needed to accont for scale separation between transport 
and droplet interactions. 

Velocity enhancement diagnosed from the model is 
non-symmetric, with larger values for sizes smaller than 
size under consideration. The asymmetry, however, may 
be related to the microphysical grid, because in E,, av¬ 
eraged within a bin parcel velocity difference is normal¬ 
ized by a difference in the terminal velocity for bins these 
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parcels belong to. As a result of the averaging and nor¬ 
malization E,, also depends on the number of bins used to 
represent domain in a radius space. Figure [TT] shows E,, 
for a REF case together with an additional 2 simulations 
using the same setup, but one with 55 bins (p=2V2), and 
other with 108 bins (p=2V4). With the increasing num¬ 
ber of bins E,. also increases, reaching the value of 900 
when the smallest bin is under consideration. A E,, values 
approach 1 independent on resolution in a radius space 
when the differences in droplet sizes is large. 

Turbulent velocity enhancement can affect signihcantly 
drizzle formation. Simulation BIN, where droplet vertical 
velocity has been set to a terminal velocity, based on the 
bin to which droplets belong produces negligible amount 
of rain water, and much smaller than for other simula¬ 
tions droplet sizes. Note, however, that even in this case 
droplets do collide forming larger ones. 

For simulations including a representation of air tur- 
bulen velocity fluctuations, drizzle forms initially prefer¬ 
entially near the cloud edge, near the entrainment eddies 
or near the cloud top, and often in the areas where is 
elevated. The edge of the cloud is a place where due to 
entrainment and mixing droplet spectrum is broader than 
in the center of the cloud and as a result coalescence be¬ 
tween droplets is more efficient. Additionally, because of 
the gradient of the velocity near the cloud edge, the dif¬ 
ferences in the relative velocity of parcels there is large 
and this also enhances probability of droplet coalescence. 
Formation of the hrst drizzle near the cloud edge it a bin 
model with the representation of the effect of turbulence 
on droplet collision rate h as been also reported recently 
by Benmoshe et al.l ( 2012 ). 

The cloud droplet spectrum is relatively broad (of the 
order of ~ \^m) from the onset of the cloud formation. 
This value is much larger than reported for parcel mod¬ 
els in the past and large enough to trigger the coalescence 
process even for the case with a very high cloud droplet 
concentration. It’s acknowledged, however, that when the 
standard deviation for each computational grid is consid¬ 
ered values smaller than 0.1 /rm are observed, indicating 
that very narrow (comparable to a parcel model) droplet 
distributions also form within computational grids. 

Coalescence does process aerosol, but the time scale of 
this process is longer than ~20 minutes of cloud evolution 
discussed here. Processing by multiple clouds is needed 
for the boundary layer aerosol to change the shape of the 


aerosol distribution. 


5. Appendix 


Consider conservation equation for a scalar for an 
anelastic approximation; 


dpC 

dt 


J— = -V(pyC) H- V(pAVC), 


( 12 ) 


where C - is concentration, p - gas density, v - gas velocity, 
K - diffusion coefficient (note that in general A is a ten¬ 
sor, however, here it is assumed that this tensor has only 
diagonal values not equal to 0; K,p=b,pW). This equation 
can be written as: 


^ = -V[pC(v + -V(pA)] + A(pAC). 
ot p 


(13) 


This equation corresponds to the following Stochas- 
tic Differentia l Equa tion in Ito sense (see for insta nce: 
Salamon et al.l ( 2006 ) or Spivakovskava et al. ( 2007 ) for 
the case when p is constant): 

dx = iv+ -X(pK))dt + y/^dW(t), (14) 

p 

where dWit) is a Wiener process. Integration of this equa¬ 
tion gives: 

1 n/S.t 

{v+ -V(pK))dt+ ^dW(t), (15) 

Jo P Jo 

or in numerical representation for i-th direction, assuming 
that V, K and p are constant during the time step; 


AA; = iViiXi) ■ 


1 


p{Xi) dxi 


(p(A,)A(A;))Af-H ilKAtG, 


(16) 

where G is a random number having normal distribution. 
It follows that effectively turbulent diffusion corresponds 
to a random walk process with mean velocity having two 
components, hrst being the velocity from the numerical 
model, the other accounting for a change in the density 
and diffusivity in a non-homogeneous medium: 


Ve - Vi(Xi) + 


1 


p(Xi) dxi 


(p{Xi)K(Xi))- 


(17) 
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and a random component: 



To find K in numerical model an assumption must be 
made about the mixing length l,K ~ l-T Typically 

it is assumed to be of the order of a grid length. Although 
in the Eulerian model this assumption is justified, because 
mixing within each model grid must be completed within 
the time-step, in a Lagrangian transport model, where the 
exact location of each parcel is known it’s not necessary 
true. For a Lagrangian model a different length scale can 
be used, for instance from model time-step and turbulent 
kinetic energy (TKE), which is a prognostic variable in a 
sub-grid scale model. The following length scale can be 
defined: li=At y/TTKEjl, which for a 2D case discussed 
in this article is equal to '^TKE. It follows that 

Vr ~ ^TKE. (19) 

This length scale is used only in a Lagrangian model to 
calculate diffusion coefhcient and next the standard devi¬ 
ation of the velocity fluctuations and mean deterministic 
velocity. 

Diagnostic from the model indicate that the velocity as¬ 
sociated with the variability in K and p is very small ~ 1 
mm/s, but this term is nevertheless kept in the model. 
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Figure 2: Upper row - cloud water mixing ratio (qc) after 8 min. (a), 15 min (b) and 17 (c) min.; lower row - rain water mixing ratio (q^) after 15 
min. (d), 17 min (e) for a REF run. Solid black line shows contour for value 0 of the vertical velocity. 
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Figure 3: Evolution in time of a Total Water Path (a) and a Rain Water Path (b) for REF - blue, BIN - red, TURB - green, REF2 - black. 
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Figure 4: Vertical profiles for a REF case of mean radius (a), standard deviation of the droplet spectrum on given level (b), (c), and cloud droplet 

concentration (d) for times 6 min. - blue, 8 min. - green, 15 min. - red. 
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Figure 5: As figure|4] but for a BIN case. 
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Figure 6: Relation between aerosol size and droplet size for grids with qc > 10 ^ - Logio(N(r,ra)) for run REF (a) and BIN (b). 


18 














dN(r)/dlog(r) [cm^ dN(r)/dlog(r) [cm 




r[nm] 



r[nm] 





r[|im] 



r[nm] 


Figure 7: Droplet spectrum after 8 min. (a), 12 min. (b), 15 min (c) and 18 min (d) for a REF - blue, BIN - red, TURB - green, REF2 - black. 
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Figure 8: As figure[7]but for aerosol spectrum. 
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Figure 9: A velocity enhancement factor (Ev.)for a different microphysical bins: bin 4 (a), bin 10 (b), bin 16 (c), bin 22 (d) for a REF - blue, BIN - 
red, TURB - green, REF2 - black. 
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Figure 11: A velocity enhancement factor (Ey) for a different microphysical grids: for radius 2.6 fim panel a, 9.1 fim - panel b, 28.8/im - panel c, 
94.3 - panel d, for simulation: REF2 - black, REF2 with 55 bins used to represent droplet domain - red, REE2 with 109 bins used to represent 

droplet domain - blue. 
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